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ABSTRACT 

In the present work we describe the formalism necessary to derive the properties of dark matter halos beyond 
two virial radius using the spherical collapse model (without shell crossing), and provide the framework for the 
theoretical prediction presented in Prada et al. (2005). We show in detail how to obtain within this model the 
probability distribution for the spherically-averaged enclosed density at any radii P(5, r). Using this probability 
distribution, we compute the most probable and mean density profiles, which turns out to differ considerably 
from each other. We also show how to obtain the typical profile, as well as the probability distribution and 
mean profile for the spherically averaged radial velocity. Three probability distributions are obtained: a first 
one is derived using a simple assumption, that is, if Q is the virial radius in Lagrangian coordinates, then the 
enclosed linear contrast Si(q) must satisfy the condition that Si(q = Q) = 5 V j r , where <5 V! > is the linear density 
contrast within the virial radius /?,,,> at the moment of virialization. Then we introduce an additional constraint 
to obtain a more accurate P(5,r) which reproduces to a higher degree of precision the distribution of the 
spherically averaged enclosed density found in the simulations. This new constraint is that, for a given q > Q, 
Si(q) < S V i r . A third probability distribution, the most accurate, is obtained imposing the strongest constraint 
that Si(q) < S V i r V q > Q, which means that there are no radii larger than /?„> where the density contrast is 
larger than that used to define the virial radius. Finally, we compare our theoretical predictions for the mean 
density and mean velocity profiles with the results found in the simulations. 

Subject headings: cosmology:theory — dark matter — large-scale structure of universe — methods: analytical 



1. INTRODUCTION 

The study of the density profile of cold dark matter halos 
beyond the virial radius is a subject of considerable relevance. 
From an observational point of view, knowledge of the shape 
of the density profile far beyond the virial radius is essen- 
tial for an appropriate interpretation of gravitational lensing 
phenomena (e.g. Smith et al. 2001; Guzik & Seljak 2002; 
Hoekstra et al. 2004; Sheldon et al. 2004), the pattern of Ly- 
man alpha absorption around virialized systems (e.g. Barkana 
2004; Bajtlik, Duncan & Ostriker 1988) as well as the mo- 
tion of satellite galaxies as a test for dark matter distribution 
at large radii (Zaritsky & White 1994, Zaritsky et al. 1997; 
Prada et al. 2003, Brainerd 2004; Conroy et al. 2004). From 
the theoretical point of view, the study of the properties of 
dark matter halos at several virial radius in cosmological sim- 
ulations provides an excellent benchmark for developing and 
testing the basic theoretical framework which will be decisive 
for a full understanding of the physical origin and formation 
of the ACDM halos. 

Understanding halo properties involves a set of theoreti- 
cal considerations. First, we have the issue of choosing the 
correct initial density profile. Also, there is the question of 
which processes are relevant to the gravitational evolution of 
the initial profile: is only the spherical collapse what matters 
or is triaxiality important? up to which radius can we use the 
standard spherical collapse model without shell crossing? are 
highly asymmetrical processes, like merging, relevant? In or- 
der to answer these questions it is very convenient to focus 
first on those properties of the halos which involve the fewest 
theoretical uncertainties. 

The dark matter density profiles at several virial radius 
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are particularly suitable to check whether the spherical col- 
lapse model can provide accurate predictions (see Prada 
et al. 2005). In fact, it has been shown that the spher- 
ical collapse model reproduces very well the relationship 
between the small values of the spherically-averaged en- 
closed density at those large distances and the radial velocity 
dLilie et al. 19861) . 
We define the spherically-averaged enclosed density as: 



p(<r) 

< Pm > 



1 + 5 



where S is the enclosed density contrast and < p„, > the av- 
erage matter density in the Universe. We can also define the 
spherically-averaged local density as: 



< Pm > 



1 + 6' 



where 8' is the density contrast in a narrow shell of radius r. 
We can then obtain the density contrast 5' from the enclosed 
density contrast 5 using the relation: 

3r z ar 

Despite to all the effort done to understand the central dense 
regions of the dark matter halos in cosmological simulations, 
not much attention has been devoted to the study of the re- 
gions beyond the formal virial radius, i.e. the radius within 
which the spherically-averaged enclosed density is equal to 
some specific value. The main goal of the work presented in 
this paper is focused on the outskirts of the dark matter halos, 
where the correct evolution of the spherically-averaged en- 
closed density profiles can essentially be obtained using the 
standard (without shell crossing) spherical collapse model. 
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This model, first developed by Gunn & Gott (1972) and Gunn 
(1977), describes the collision-less collapse of a spherical per- 
turbation in an expanding background. They introduced for 
the first time the cosmological expansion and the role of adia- 
batic invariance in the formation of individual objects. Later, 
Fillmore & Goldreich (1984) found analytical predictions for 
the density of collapsed objects seeded by scale-free primor- 
dial perturbations in a fiat universe. Hoffman & Shaham 
(1985) generalized these solutions to realistic initial condi- 
tions in fiat and open Friedmann models. Some studies have 
been done to include more realistic dynamics of the growth 
process (e.g. Padmanabhan 1996; Avila-Reese, Firmani & 
Hernandez 1998; Lokas 2000; Subramanian, Cen & Ostriker 
2000). 

There are plenty of works in the literature using the spheri- 
cal collapse model to predict the density profiles of dark mat- 
ter halos mainly focused on explaining their central regions. 
For example Bertschinger (1985) used the spherical collapse 
with shell crossing to obtain the density profiles resulting 
from initial power law density profiles. Lokas & Hoffman 
(2000) considered more general initial profiles. The effect 
of non-radial motions has also been widely treated (see Ry- 
den & Gunn 1987; Gurevich & Zybin 1988; Avila-Reese et 
al. 1998; White & Zaritsky 1992; Sikivie et al. 1997; Nusser 
2001; Hiotelis 2001). Some of these authors have used ar- 
bitrary initial profiles, while others have assumed the mean 
initial profile around density maxima (Bardeen et al. 1986, 
BBKS). In all these works angular momentum is introduced 
by hand, although more recently it has been done in a more 
natural way (Nusser 2001; Ascasibar et al. 2004). Concern- 
ing to the outer parts of the dark matter halos, only Barkana 
(2004) has adopted an appropriate initial profile, but only for 
a restricted type of density profile (the typical profile). The 
more recent work by Prada et al. (2005) have obtained pre- 
dictions for the mean and most probable density profiles and 
have provided a detailed comparison with cosmological sim- 
ulations. 

A proper understanding of the physics of dark matter halos 
involves predicting correctly not only the mean halo density 
profile for any given mass but also the whole probability dis- 
tribution for the enclosed density contrast at any given radii, 
P(5, r). A first attempt to determine it can be found in Prada 
et al. (2005), where it has been shown to be generally in good 
agreement with the cosmological simulations. Nevertheless, 
this probability distribution shows, at any radius, a longer tail 
for large values of <5, as compared to that from simulations at 
any radius. In this paper we present a more accurate predic- 
tion for the probability distribution P(6, r) that constitutes the 
main new result of this work. We also give in detail the the- 
oretical background of the predictions presented in Prada et 
al. (2005). The agreement of our new predictions with the 
simulations is excellent even in the tail of the distribution. 
Furthermore, we also compute the radial velocity probability 
distribution and the mean radial velocity profile. 

The work is organized as follows. In section 2 we present 
our theoretical framework and obtain the typical density pro- 
file of dark matter halos. In section 3 we show in detail how to 
obtain the probability distribution for the spherically-averaged 
enclosed density contrast at a given radii, P(8, r), presented in 
Prada et al. (2005). The most probable and mean profiles are 
derived. In section 4, we compute the probability distribution 
and mean profile for the spherically averaged radial velocity. 
In section 5, we obtain more accurate probability distributions 
than that used in the previous sections, and compare again 



with that found in the simulations. Final remarks are given in 
section 6. 

2. THE TYPICAL DENSITY PROFILE OF DARK MATTER HALOS 

The present spherically-averaged enclosed density profiles 
attains a density contrast value of A V! > at certain radius, the so 
called virial radius. At larger radii the density contrast must 
be, by definition, smaller than A,,,,-, otherwise the virial radius 
would be larger than its nominal value. 

We shall now make some comments on the values of A vir 
and 5 V i r that we use: although at several virial radius the 
spherically-averaged enclosed linear and actual densities are 
related by the spherical collapse model, the same does not 
apply within the virial radius. The spherically-averaged en- 
closed density contrast within one virial radius, A V! >, and the 
corresponding enclosed linear density contrast, <5 vir , are not re- 
lated as homologous quantities at larger radii, because at one 
virial radius shell crossing has already becomes important. 
Consequently, the value of <5,,,> corresponding to A,,,, =340 (the 
value we adopted to define the virial radius) is somewhat un- 
certain. As a result of work still in progress we will be able 
to provide the precise values for <5 v ; r and determine its possi- 
ble small dependence on mass. Here we use for all masses 
(5 v ,v=1.9, a value that leads to very good results and that may 
be inferred from the fact that whe n A v ,>=180, 6 V j r see ms to be 
close to 1.68 for all cosmologies (Jenki nTet al. 200lh . 

It must be noted that for all our predictions it is irrelevant 
whether the value of A V! > that we use actually corresponds to 
the virial density contrast or not. By virial density contrast is 
usually meant the enclosed density contrast within the largest 
radii so that we have statistical equilibrium. The precise value 
of A vl > corresponding to this definition is still problematic but, 
to our purposes, it can be chosen freely to define a conven- 
tional "virial radius". Since we have used numerical simula- 
tions with A V( > equal to 340, we will take the same value by 
default in all the calculations. 

Let Si{q,i) be a realization of the spherically-averaged ini- 
tial enclosed density profile around a protohalo (with a given 
present virial radius, R) linearly extrapolated to the present, 
where q is the lagrangian distance from the center of the 
halo to the given point and i is an index running over real- 
izations. Any realization of the initial profile may be trans- 
formed using the standard spherical collapse model (without 
shell-crossing). We can use the relationship between the lin- 
ear value of the density contrast within a sphere, Si, and the 
actual density contrast within that sphere, S, for the concor- 
dant cosmology (Sheth and Tormen (2002)): 



Si(S)-- 



1.676 
1.68647 



1.68647- 



1.35 



1.12431 0.78785 



(l + <5) 2 /3 (1 + <5)V2 (1 + 5)0.58661 

(1) " 

or, rather its inverse function 5(5i) (Patiri et al.(2004) ex- 
pression (4)): 



( 5((5,) = 0.993[(l-0.607((5,-6.5 x 10" 3 (1 -<9(<5,)+ 



+ 0(4 -1.55))$))- 
being 9 the step function: 



1] 



(2) 



9{X): 



1 ifx>0 
ifx<0 



Transforming for every shell 6/ and q into S and r (the Eu- 
lerian radius of the shell) we may obtain, in parametric form, 
the initial profile Si(q,i) spherically evolved, 5(r, i). 
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8{r,i) = 8{8 l {q,i))\ r=[l + S{S,{q,iW' /3 q (3) 

This two equations gives 8(r, i) implicitly with q as param- 
eter. 

We must now eliminate all linear profiles leading to present 
profiles which attains an enclosed density contrast larger than 
(or equal to) A V! > at a radii larger than the nominal virial ra- 
dius. The ensemble of remaining halos allow us to obtain the 
predictions of the standard spherical collapse model for any 
statistics. For example, we shall obtain the predictions for the 
mean profile: 

8(r)=<8(rJ) >j 

That is, the average over all remaining halos (J runs over 
these halos). We shall also consider the most probable profile, 
S p (r), that is, the profile that associates with every value of r 
the 8 value with the largest probability density. 

The procedure described above to obtain the spherical 
model predictions (i.e. realizations of the linear profile 
evolved with expression (0) serve mostly to the purpose of 
clarifying the meaning of those predictions and we only use 
it as a test to the analytical expressions. In practise, we shall 
use another procedure to obtain directly S(r), 8 p (r) and, in 
fact, the whole probability distribution for the value of 6 at a 
given value of r, P(8, r). 

Before dealing with the detailed predictions just mentioned 
we consider a simpler prediction which shall help clarifying 
the rest o f the work, an d which is related to previous ap- 
proaches dBarkana 20041) . 

Consider the ensemble of all halos, 8(r,j), attaining a 8 
value equal to A V! > at a virial radius and smaller values for 
larger radii. If we transform back these profile to their linear 
counterpart, we obtain the ensemble Si(q,j). Let us now take 
the average over this ensemble (now for a fixed value of q): 

8(q)=<S(qJ) >j 

Evolving this profile by means of the spherical collapse 
model we obtain a profile which we call typical profile and 
represent by 8,{r), that is, this profile is simply the mean pro- 
file in the initial conditions spherically evolved. 

Note that the typical density profile is defined because of its 
simplicity and not because it constitutes a prediction for any 
specific statistics of the actual halos. Nevertheless, it should 
not be very different from the most probable profile. Later we 
will study how these two profiles differ from each other. This 
profile definition is the same as that used by Barkana (2004), 
but we use a different approximation to derive it. 

To obtain the mean linear density profile 6i{q) subject to 
the condition that in the present enclosed density profile the 
virialization density contrast, A,,, r , is not attained beyond the 
virial radius, Barkana used barrier penetration results. He ob- 
tained the probability distribution for <5/(q) given the condi- 
tions: 

Si(Q) = 5 vir ; 5i(q) < S vir Vq>Q 

Q = R vir (l+A rir )^ 3 

where R v , r is the virial radius (Q the corresponding La- 
grangian radius), and <J V! > is the linear counterpart of A w >. 
Then, by averaging over all 8/ values smaller than 5 V! >, he ob- 
tained the mean linear profile, Si(q), with those conditions. 
However, he had to make some simplifications, the most rel- 
evant one is that he used a sharp filter in k-space, rather than 



the top hat filter in ordinary space which is the natural one in 
this context. 

In our approach we first obtain the probability for Si{q) only 
with the condition S/(Q) = S V i r : 



exp 

P(5 u q) = P(S,(q)\5,(Q) = S vir ) = 

where 

g(q) = er 2 - -^2— 



2tt g? 



(4) 



<ri = *(G); 
1 



2tt 2 



CTi2 = cr n (q) = 



2n 2 J 



02 = cr(q) 



8 k | 2 W}(xk) k 2 dk 



8 k | 2 W T (qk) W T (Qk) k 2 dk 



W T (x) = 



3(sin x —x cos x ) 



where |<5a | 2 stands for the power spectra of the density fluc- 
tuations linearly extrapolated to the present. 

It is convenient to use a simple and accurate approximation 
for aniq)- 



crn(q) 



-b(Q)((p 2 -i) 



(o-(0) 2 "" (5) 

where b(Q) is a coefficient depending on the the size of the 
halo, Q: 



2 dmx 



x=Q 



If no restriction other than (5/(0 = <5wr were imposed on 
S/(q) the mean linear profile would be: 



Si(q)- 



P(Sl(q)/Si(Q) = S vir ) 6,{q) d(S,(q)) = 5 vir ^^ = S (q) 



However, we are interested on the mean profiles satisfying 
also Si(q) < S v i r for q > Q. So, we must use as mean profile: 



S,(q) = 



jZmWm = 5nr) S,(q) d{8,{q)) 



(6) 



= 8 (q)- 



a(q)e 2 



1 - ^er 



.f r ( S„-S {q) \ 



<?(q) = (g(q)) 



1/2 



For arbitrarily massive halos a(Q) << 1. So, since for the 
relevant q values (q > Q) W{q) < a{Q), 8i(q) is simply given by 
8o(q). However, for halos with a(Q) > 1, 8i(q) is substantially 
steeper than 8o(q), resulting in steeper present density pro- 
files for smaller ma sses. This resu lt has previously been ad- 
vanced by Barkana (Barkana 2004) and have been co nfirmed 
by means of numerical simulations dPrada et al. 20051) . 
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It must be noted that the profile given by expression (0 is 
not exactly the mean linear profile implicit in the definition of 
typical profile. Note that at each value of q, Si is constrained 
to lay below <j v ,,> but the probability distribution upon which 
this constraint is imposed does not account for the fact that 
the profile lies below (5,.,y at any other value of q larger than 
Q. This is, however, a good approximation, because the most 
relevant part of the present density profile (to k 10 virial ra- 
dius) corresponds to a narrow region in Lagrangian coordi- 
nates (w 1 .50. The value of 5/ for q between Q and 1 .5Q are 
strongly correlated. So, if we impose the condition Si < <5,,,y 
at, for example, q = 1 .25Q, the probability that the same con- 
dition holds at any q is close to one. 

Once we have 5[(q) all we need to do is to evolve it with the 
spherical collapse model. Using equations we may write 
(Patiri et al.(2004) expression (20)): 



S(r) = S(5 l (q)); q = r[l + 8(r)]i 



(7) 



where the right hand side of this equation is simply the 
function defined in expression (0 evaluated at Si(q) (given 
by expression (0). For each value of r we must solve this 
equation for the variable S(r). Applying the function defined 
in expression (0, which is the inverse of that defined in (0, 
to both sides of this equation we have: 



5 l (5(r)) = 5,(q); q = r(l + S(r)y 



(8) 



where the left hand side is expression (0 evaluated at S{r). 
This equation is usually simpler to solve than equation (0 and 
is the one we used in Prada et al.(2005). It must be noted, 
however, that, if one intends to generate the whole profile 
rather than to obtain S(r) for some specific value of r, it is not 
necessary to solve equation (0, since one may obtain all cou- 
ples of values of S{r), r using (0 and running over all values 
of q. The profile obtained in this way is the typical enclosed 
density contrast profile, 5,(r). To obtain the typical density 
contrast profile, 5'(r), we may use the relationship ( TToT l. given 
at the end of next section. 

3. THE PROBABILITY DISTRIBUTION, P(8,R). MOST PROBABLE 
AND MEAN PROFILES 

At a given value of r, S takes different values, S(r,j), over 
the assemble of halos. The question now is which is the prob- 
ability distribution of S over this ensemble, P(S,r). As we 
saw in the previous section, this can be done, in principle, by 
making realizations of the initial profile, S[(q,i), and evolv- 
ing them accordingly with equations 0. Let's assume, as a 
first approximation, that the realizations of the initial profile 
may be carried out by generating for each value of q a value 
of Si accordingly with distribution (0. That is, we assume 
that the distribution of Si is only conditioned by the fact that 
S/(Q) = 5 V j r . We shall latter consider initial profiles with an ad- 
ditional constraint. With these realizations we can elaborate 
for each value of r a histogram for P(S,r). There is, how- 
ever, a direct analytical procedure to obtain P(S,r) from the 
probability distribution for 5/ (expression (0 in the present 
approximation). 

S is a unique function of 5/ (expression (0). So, one may 
think that P(<5,r) can be obtained from simply through the 
change of variable S = 5(5[). However, expressions (0 show 
that in transforming the initial profile not only is 5/ trans- 
formed into 5, but also q is transformed into r. Now, expres- 
sion (0 gives the distribution of 5/ at a fixed q value. But what 
we want to obtain is the distribution for 5 at a fixed r value. 



So, since the relationship between q and r depends on 5 (or 5i) 
itself, it is clear that the derivation of P(5,i) from expression 
(0 (i.e. from P(<5/,q)) can not be as simple as described above. 

Fortunately, there is a simple expression relating both prob- 
ability distributions, which is valid as long as shell-crossing is 
not important: 



P(S,r) = - 



dd 



f_ZP(5 h q) dSi 



(9) 



q = r (l + <5) 3 



where Si (5) is given by expression and P(S/,q) is the 
linear profile. Note that S enters not only in the integration 
limit but also in the integrand through q. The derivation of this 
relationship is given in Patiri et al.(2004) appendix B. In that 
work this relationship was derived in regard with void density 
profiles, so it had an slightly different form. In appendix B we 
give the derivation corresponding to the present case. 

Using expression (0 for P(<5/,q) we find: 

&) P ( 6 ' >d d6 ' 1 erfc(F(x = 5,(5))) - erfc(F(x = S vir )) 



f_-p(S h q)dS, 



1 - \erfc(F(x = S vir )) 



(10) 



F( X ): 



<j(Q) v " 



q = (r(l + 6))i 

with (Jniq), cr(Q), g(q) as defined in (0. 

For the purposes of this section we may neglect the term 
erfc(F(x = <5,,,>)). The full expression shall be used in section 
5 along with a more refined one. We then have: 

P(S, r) = 4= e- F2(x = s >< 5 » ^-F(x = 5,(5)) (11) 
y/ir dd 

By construction, 5 at r must be smaller than A V! > and larger 
than certain value, <5, m „(r): 



<Wr) = 341 (R vir /r) 3 -l 



(12) 



This minimum value corresponds to a situation where there 
is no matter in between R V j r and r. So, for 5 values outside the 
interval (5 mi „(r), A,,,>), P(5, r) is zero. 

We may obtain an analytical expression for P(<5,r) using ap- 
proximation @ for and the following approximation for 
<r(q) (which enters g(q), defined below expression (0): 

(7(4)^(1.65 \Qr 2 +0.\05{qhMpc~ x )T^ for q < 3 IT x Mpc 

As mentioned before, expression (0 is valid as long as 
shell-crossing is not important. In Prada et al.(2005) we 
have found by mean of comparison with numerical simula- 
tions that, beyond three virial radius, the relevance of shell- 
crossing diminishes quickly. This relevance can be estimated 
a priori (i.e. without comparison with simulations) obtain- 
ing P(<5,r) directly through realizations of the initial profile ac- 
cordingly with expression (0 and evolving them accordingly 
with equations (0. If shell-crossing were irrelevant, the P(<5,r) 
obtained in this way should be equal to that given by expres- 
sion ( fTTT ). The presence of certain amount of shell-crossing 
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will cause the P(<5,r) obtained with realization to have a some- 
what smaller maxima and a more extended tail than that given 
by ([TlT l (note that it is this expression that corresponds to re- 
alizations elaborated with expression (O. In this case none of 
the procedures gives the correct P(<S,r) because both assume 
that 8 is related to 81 by means of expression (Q~|), which is 
inconsistent when shell-crossing is important. However, the 
difference between the results obtained with both procedures 
is of the same order of the difference between any of them and 
the profile obtained with a proper treatment of shell-crossing. 
Note that even this last P(<S,r) is not the real one, since, as we 
said before, we are generating the initial profile using only a 
two point distribution (expression 

In figure [TJ we compare the P(<5,r) obtained by the two pro- 
cedures mentioned above for several values of r expressed in 
unit of the virial radius (that we denote by s). We also show 
the corresponding histograms obtained from the numerical 
simulations described in Prada et al.(2005), which were done 
using the Adaptive Refinement Tree (ART) code (Kravtsov et 
al. 1997) for the standard ACDM cosmological model with 
= 0.3, J7a = 0.7, h = 0.7, and cover a wide range of scales 
with different mass and force resolutions (see Prada et al. 
2005 for a detailed description). In particular, the histograms 
in figure 1 were obtained for a total of 654 halos in a mass 
range of 3 ± 1 ■ I0 l2 h~ l M selected without any kind of iso- 
lation criteria. 

We can see that, for s = 2.5, where shell-crossing is already 
important, there is substantial difference between the results 
of both procedures. A considerable amount of probability is 
transfered from the most probable value to much larger values 
(5 ps 100) causing the distribution obtained through the real- 
izations to be bimodal. For s = 3.5 there is still a small amount 
of shell-crossing causing the maxima obtained with both pro- 
cedures to differ by roughly a 20%. For larger values of s this 
difference steadily diminishes. 



o.ic 



0.14 






FIG. 1. — P(S,r) as given by expression (TTJ compared with the corre- 
sponding histogram obtained through realizations (histogram with small er- 
ror bars) and that found in the simulations (histogram with large error bars) 
for four values of ,v (= r/R r j r ) and a mass of 3 X 10 12 /r ! M Q (A vi , = 340 
and5 v „ = 1.9). 



with the numerical simulations reduces. In fact, even for 
s = 3.5 the relative values of P(8,i) at different values of 8 
to the left of the maxima are very well given by (fTTT >. The 
difference in the absolute values with respect to those in the 
histogram is due to the normalization. For 5 values to the right 
of the maxima, expression (fTTT i gives a considerably more ex- 
tended tail than in the actual distribution. Therefore the nor- 
malization constant is larger in the latter case. 

Note that, although the maxima of expression (flTT i and that 
of the histogram corresponding to the numerical simulations 
approach as s increases, due to the increasing irrelevance of 
the tail, this tail is still substantially more extended even for 
5 = 6. Now, since the relevance of shell-crossing is small for 
s larger than 3, the most likely explanation for the excess in 
the tail given by expression ( fTTT i lies on the fact that we are 
using expression (01 for P(<5/,q). In the last section we shall 
consider a better P(S,i) and discuss the resulting improvement 
of the behaviour of its tail. 

Using P(<5,r) (as given by (fTTt ) we may immediately obtain 
the most probable and the mean profiles. For the first one we 
have: 



d p (r) = 5 max ; -yxP(S,r) 
do 

And for the mean profile, in principle: 



8=8,,, 



P(S,r)5 dS 



(13) 



(14) 



However, due to the fact that the mean is rather sensitive to 
the form of the tail, we must artificially cut off the tail. From 
the simulations we know that the real tail practically ends at 
S ~ So(r) with So given by: 

P(8 maX7 r) 



P(5o,r) = 



25 



So, instead of (|T4]> we use (Prada et al.2005): 



S(r). 



J^ } P(S,r)SdS 
J^ r) P(5,r)dS 



(15) 



In table 1 we give the values of So and S max for several values 
ofi(= r/R vir ). 

So far we have considered the probability distribution and 
profiles of the spherically-averaged enclosed density contrast, 
S. We shall now consider the local density contrast S\ that 
is the density contrast within a narrow shell of radius r. In 
this case we can not obtain P((5',r) by mean of a simple rela- 
tionship like expression ([9j. However, the mean 8' profile can 
be obtained from the mean 8 profile. To this end, note that 
for each actual density profile, that is, for each spherically 
evolved realization of the linear profile, S(r,j), the following 
relationship holds: 



S'(rJ) = ^2 ^-(r i 8(r,j)) 
br L dr 



(16) 



which follows immediately from the definitions of S, 8'. 
The mean 8 profile does not correspond to any actual density 
profile. However, being the mean a linear operation, the same 
relationship holds for the mean profiles: 



Furthermore, as s increases the difference between the dis- 
tribution given by expression ( fTTT i and the histogram obtained 



1 d 



8 (r) =< S'(rJ) >j=< — —(r 3 S(rJ)) >j 
3r z dr 



(17) 
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TABLE 1 
Artificial cut-off for P(5, r) 



s 


&max 


<5o 


1.5 


115.2 


495.4 


2.5 


28.6 


135.3 


3.5 


11.3 


67.9 


4.5 


5.8 


42.2 


5.5 


3.3 


29.3 


6.5 


2.0 


19.6 


7.5 


1.3 


12.8 


8.5 


0.82 


8.1 



The most probable profile is constructed with the ensemble 
of profiles S(r,j) by means of a non-linear operation: that of 
choosing, for any value of r, the member of the ensemble with 
the largest probability density for S(r). So, relationship ( fT6b 
is not valid between S' p and S p , because the value of S p at 
different values of r may correspond to different members of 
the ensemble. However, in Prada et al.(2005) we have used 
expression ( fTSl l to obtain S' p (i) from S p (r) and found results in 
good agreement with the simulations, but this is a posteriori 

agreement: unlike the prediction for 5 (r), the value of S' p (i) 
obtained using expression ( TT6T > is not a proper prediction. 

4. RADIAL VELOCITY PROFILE OF DARK MATTER HALOS 

In the spherical collapse model, the radial velocity at dis- 
tance r with respect to the center of the spherical cloud is a 
unique function of the spherically-averaged enclosed density 
contrast, 6{r). An exact analytical expression may be given 
for this relationship using 5i(5) (expression (HJ). Mass conser- 
vation within a shell with initial Lagrangian radius q, implies 
that, at any time, the following relationship must hold: 



K0 

a(t) 



(18) 



where a(t) is the scale factor of the universe (normalized to 
1 at present) so that the left hand side is the comoving radii 
of that shell at time t . Deriving this equation with respect to 
time we find: 



a(t) a\t) ' ' a(t) 3 H 

where H{t ) is the Hubble constant. Using now 5{5j) we may 
express 5 in the form: 



• D(t) a(t) dD{a(t)) 

Oi = 0; = H Oi (20) 

D(t) D(a(t)) da(t) 



where D{t) is the growing mode of linear density fluctu- 
ations. Now, since <5(<5/) is the inverse function of Si(5) we 
have: 

dS(Si) fdSiWY 1 



dSi 



dS 



Using this in (|20| >, inserting (|20T > in ( fl9] > and using again 
( fT~8T > we find: 



1 a(t)H dD{a(t)) r 5,(5) 
3D(a(t)) da(t) {1 + S) 



dS,(S) 
dS 



Vr = f = H(t)r-- 1 ^—^^— R ^ = rf(S) 

(21) 

For the concordant cosmology we have at present H = 
72 km/s/Mpc, o=l, ^f/D{a) = 0.51. Writing (ED in the 
form V r = rf(8) we may use expression ( fTTT i to obtain the 
probability distribution, P(V r , r), for V r at radii r. 



P(V r ,r) = P[S = r 



(22) 



where f~ l is the inverse function of / and where P(S,r) is 
given by expression ( fTTT i. 

This distribution should not be confused with the distribu- 
tion for the radial velocity of dark matter particles at a given 
value r. V r is the mean radial velocity of all particles in a given 
narrow shell with radius r. So, for a given halo and a given 
value of r, V r takes a unique value. Expression ( l22l gives the 
distribution of this value over the ensemble of halos. 

The mean V r profile is given by: 



f 



So(r) 



P(6,r)rf(6) dS 



S!^P{5,r)d5 



(23) 



with P(8, r) given by (fTTT i. 6q as given in Table 1, f(6) given 
by (ED an d Smin given by (fT2l . 

In figure E] we show the predictions given by d23l i for the 
mean radial velocity profile, and compare it with the profile: 



V r = rf(6(r)) 



(24) 



That is, for each value of r, this expression gives V r corre- 
sponding to the mean 6 at that r through Eq.(t2TT>. We find that 
Eq.d24b is a good approximation to Eq.([23l. 

Note that, although we have used the probability distri- 
bution given by Eq. dTTl) . these expressions are valid for any 
P(S,r). 



200 - 




FIG. 2. — Radial velocit y pro file as given by expression (23) (filled line) 
and as given by expression (24} (dashed line). Both correspond to a mass of 
3 X 10 12 /r'M 
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5. IMPROVING P(6, R) 

We have seen that for s (= r/R V j r ) larger than 3 the realiza- 
tion of initial profiles and their subsequent evolution lead to 
values of P(S,r) very close to those obtained using expres- 
sion ( fTTT i. This means that shell crossing is not important at 
these radii. So, the difference between the actual value of 
P(5, r) (the histogram obtained from the simulations) and that 
given by expression (fTTT i lies on the fact that this expression is 
only an approximation, or in the possible relevance of triaxi- 
ality, non-radial motions and pressure (velocity dispersions). 
To determine the amount of the discrepancy due to inaccura- 
cies in the statistical description of the initial conditions as op- 
posed to the discrepancy due to neglected dynamical factors, 
we must consider better approximation than that provided by 
expression (fTTT i. We shall use now two improved approxima- 
tions. The first one is that obtained by using expression ( TTOb 
in expression (O, which we represent by P\(S,r). Note that, 
to derive this approximation we have used the fact that at a 
given q larger than Q, Si must be smaller than <5 V! >, but we 
have used expression (|4|i for P(5i,q). This distribution does 
not account for the fact that for all values of q larger than Q, 
Si < 8 V j r . Consequently, P(Si,q) falls off with Si more slowly 
than it should and the same applies to P(S,r). In the second 
approximation, which we represent by Pi{8, r), we consider a 
P(5i , q) which takes into account this additional constraint. By 
comparing the predictions obtained with PiiS, r), P\(S, r) with 
those obtained with P{5,r) (expression (fTTt ) we may check 
whether P\ (S, r) is accurate enough, so that the remaining dis- 
crepancy of the prediction with the results shown in the sim- 
ulations may be ascribed to unaccounted dynamical factors. 

We have already pointed out that the most relevant part of 
the present density profile (up to ~ 8 virial radii) comes from 
a narrow band in Lagrangian coordinates (from Q to ~ 1 .50 
so that the values of Si within this band are strongly corre- 
lated. This is the reason why using expression <j4j> for P(Si,q) 
is a good approximation, because if <5; at q (within the men- 
tioned band) lies below <5 V! > the same will probably be valid 
also at all q' larger than Q. We go now a step further and ex- 
plicitly demand that Si lies below 5,,,> at all Lagrangian radii 
in between Q and q. Imposing the same condition at q' larger 
than q is unnecessary since it will lead to a negligible change 
for P(5i,q), because if Si lies below 5,,,> at q the same will 
almost certainly occur at larger distances. 

Given the strong correlation between the values of Si in be- 
tween Q and q (for the relevant g's), the condition that for 
all values of q larger than Q Si < <5 v /r> is almost equiva- 
lent to demanding that Si lies below <J VJ > at the middle point 
q' = j(q + Q). We might have chosen any other point in be- 
tween and searched for the point imposing the strongest con- 
straint, since the real constraint must be stronger (i.e. the tail 
of the distribution falls off more steeply) than that imposed by 
the point leading to the strongest constraint. We have chosen 
the middle point because it seems a priori a good choice. 

Here we discuss the main lines of the derivation and the re- 
sult for P{8 1, q) within this new condition, which we represent 
by P2(Si,q). We give in appendix [Al the details and the full 
descriptions of the expressions involved. 

To obtain P2(Si,q) we must first obtain the joint probability 
distribution for the value of Si at Q, the middle point, and at q. 
We represent by X\, X2, X3 respectively the value of Si at these 
three points. 

The joint distribution for these three variables is a Gaussian 
trivariate distribution P(Jti,JC2>*3)> which can be obtained for 









S = 3.5 






*,-A..-* 





40 
S 




FIG. 3. — Probability distribution for 8 at 3.5 and 6 virial radius for 3 X 
10 12 h" 1 At m . T he filled curve corresponds to the approximation given by 
expression (26}, the dashed cu rve to that given by full expression (9), and 
the dotted curve to expression jilt . The histogram corresponds to the same 
simulations as in figff] 




FIG. 4. — Mean 8 profile for 3 X 10 12 h M@ using the probability distri- 
bution given by expression jilt (da shed line), by full expression (5) (3dots- 
dashed line) and by expression (26) (dot-dashed line). Mean 8 obtained from 
simulations is given for comparison (filled line). In all the cases, a maximum 
value of S = 10 was used. 

a given power spectra. With this distribution we may imme- 
diately obtain the distribution of X3 conditioned to x\ = £ V! >, 
X2 < S V i r , namely: P2(Si,q). So, we have: 



ft(«,9) = 



Jl™ P ( X2 \ Xl ^ P ( X3 \ x ^ x 2) dx 2 



(25) 



where P(S v i r) X2) is the joint probability distribution for x\, 
X2 with X\ = S V j r . P2(Si,q) is the probability distribution of S/ 
at q conditioned to <5/(0 = <5,,,> and Si(q') < S vir . We have not 
imposed yet the almost redundant condition (since P2(S/,q) 
is very small for Si > <5 V! >) that Si at q must be smaller than 
S V i r . The probability distribution for 8/ with this condition, 
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FIG. 5 . — Peculiar infall velocity profil e fo r 3 X 10 12 H~ l M@ using the prob- 
ability distribution given by expression 111} (dashed line), by full expression 
{§) (3dots-dashed line) and by expression )26t (dot-dashed line). Peculiar in- 
fall velocity profile obtained from simulations is given for comparison (filled 
line). In all the cases, a maximum value of 5 = 70 was used. The straight line 
corresponds to the Hubble Flow. 



FIG . 7 . — Peculiar infall velocity profil e fo r 3 X 10 12 h~ l Mq using the prob- 
ability distribution given by expression i 1 1 1 (dashed line), by full expression 
(5) (3dots-dashed line) and by expression i26\ (dot-dashed line). Peculiar in- 
fall velocity profile obtained from simulations is given for comparison (filled 
line). For all radius, the average was calculated excluding the 20% of the 
halos with the largest <5 values. The straight line corresponds to the Hubble 
Flow. 
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FIG. 6.— Mean <5 profile for 3 X 10 12 /t'Mq using the probability distri- 
bution given by expression lilt (da shed line), by full expression (9) (3dots- 
dashed line) and by expression i26\ (dot-dashed line). Mean <5 obtained from 
simulations is given for comparison (filled line). For all radius, the average 
was calculated excluding the 20% of the halos with the largest 8 values. 

Pi(Si,q), is simply obtained through normalization (see ( |A5t ) 
The distribution for i5 at fixed r within this formalism, 
which we represent by Pi(S, r), may be obtained from P2(Si,r) 
by means of expression (0, which is valid when shell- 
crossing is not important (s > 3): 

P 2 (S,r) = -— / P 2 (6 l ,q)dS l (26) 

dO J 5,(5) 
q = r(l + 5)^ 

In figure[3]we show the probability distributions for a mass 
of 3 x 10 12 /i _1 M Q at 3.5 and 6.0 virial radius. As expected, 



P2(6,r) falls off much faster than P{S,r) being in excellent 
agreement with the simulations. The tail of P\(S,r) falls off 
sufficiently fast to give sensible results for the density and ve- 
locity profiles averaged over all possible halos (5 between <5 m ,„ 
and A v! >). However, we know that for S > 70 the standard 
spherical collapse model is not a good approximation, so we 
can not learn much by comparing simulations with predictions 
for averages over all halos. It is more instructive to compare 
the predictions for the averages corresponding to 5 values be- 
tween 6,„i„ (expression ( fT2b ) and 70 with those found in the 
simulations for the same range of S values. 

In figureHwe show the predictions given by P(S, r), P\ (5,r), 
P2(S, r) for the 5 profile averaged between <5, m „ and 70 and 
the averaged found in the simulations. In figure [5] the cor- 
responding results for the peculiar infall velocity (the Hub- 
ble Flow minus the radial velocity) are presented. One of the 
things we learn from these figures is that P2(S,r) must differ 
very little from the exact distribution, since it may be shown 
on general grounds that the difference between the predic- 
tions given by these two distributions must be much smaller 
than the difference between the predictions given by P2(8,r) 
and Pi(S,r). Another thing that we learn is that the spheri- 
cal collapse model can not be a good approximation for all 
halos included in the averaging process. The reason being 
that only the presence of some considerably triaxial haloes 
(i.e. haloes such that the density contours of their outskirts 
are very triaxial) within the averaging assemble can explain 
the fact that the results found in the simulations are somewhat 
larger than those obtained using P2(S, r). If the dynamical fac- 
tor that we have neglected were irrelevant the latter results 
must be slightly above the former. Non-radial motion and ve- 
locity dispersions preserve the order, only increasing the dif- 
ference. The reason being that when these factors are taken 
into account, while preserving spherical symmetry, Si(S) is 
steeper (S smaller for given <5/) than in the standard model (ex- 
pression (HJ), resulting, through expression ( 1261 ) and (I A81 > in 
steeper 5(r). So, only the relevance of triaxiality can explain 
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the results shown in figure |4] The effect is not large in this 
figure, since only a small fraction of the halos are affected, 
but it is quite meaningful. It is this very effect that is causing 
that the peculiar infall velocity profile predicted by means of 
P2(S, r) does not agree with that found in the simulations even 
at r/R vir > 7. 

Knowing that for some halos the spherical collapse model 
can not be a good approximation, the relevant question now 
is to determine precisely how good is it for most halos. To 
this end, at each radial bin we search for the value, S\ (r), such 
that the upper cumulative probability, as given by P2(S,r), is 
0.2, and eliminate, both in the simulations and in the predic- 
tions, halos with larger S values at that bin. By doing so, we 
eliminate many halos which simply have flatter profiles than 
average, but are otherwise sufficiently spherical for the model 
to apply. But, we are sure of having eliminated all halos with 
highly triaxial outskirts, most of them corresponding to situ- 
ations where a couple of halos lie within a few virial radius 
from each other, so that each of them will show up in the out- 
skirt of the other halo; a situation that, by no mean, can be 
described by the spherical model. The resulting profiles, as 
given by the different approximations to P{S, r) and the simu- 
lations are given in figure [6] The difference between the vari- 
ous approximations is now smaller, since the main difference 
occurs at the far tail (of P(S,r)), which has now been elimi- 
nated, but this difference increases for smaller masses. The 
profile obtained from the simulations lies now slightly below 
the best prediction, as it should be if the spherical model is 
a good approximation. In figure [7] we show the correspond- 
ing results for the peculiar infall velocity. It is apparent that 
the prediction agrees well with the simulations at large radii 
but as we go below ~ 5 virial radius there is an increasing 
discrepancy. This must be due to the unaccounted dynamical 
factors. 

Figure[8]give us another clue as to the relative effect of these 
factors. Triaxiality causes the dispersion of 5 at any given 
radii to increase, since, for a given value of 6/ triaxial evolu- 
tion gives a distribution of S values with a finite dispersion, 
while the spherical model gives a single value. So, the fact 
that the dispersion found in the simulations lies somewhat be- 
low the predictions, imply that triaxiality can not be dominant 
amongst unaccounted factors. 

In summary, the distribution P2(8, r) is very close to the ex- 
act. To most purposes the simpler distribution P\(S, r) may 
be used, their difference being small, although it increases for 
smaller values. The density profiles obtained with P2(S,r) are 
in excellent agreement with those found in simulations be- 
yond ~ 3 virial radius. The density dispersion and the ra- 
dial velocity show some discrepancy below ~ 5 virial radius, 
which clearly indicate the relevance of unaccounted dynami- 
cal factors, specially velocity dispersions. 

6. FINAL REMARKS 

The spherical collapse model describes very well the prop- 
erties of dark matter halos beyond three virial radius. This 
could seem surprising given the fact that the density contours 
around halos may be considerably aspherical and the presence 
of tidal fields. 

Nevertheless, the assumption that the spherically averaged 
density profiles at several virial radius evolves according with 



the spherical collapse model is not just a simplification in- 
troduced to make the problem tractable. The mean evolu- 
tion is given by the spherical collapse model with some dis- 
persion (for a given S{) due to triaxiality. This have been 




FIG. 8. — a corresponding to the mean 8 profiles (excluding at each radius 
the upper 20% 5 values) obtained using different probability distributions, as 
shown in fig|6](same type of lines used here); a of the mean <5 obtained from 
simulations is also given for comparison (filled line). 



shown by mean of simulations (Lilje et al. 1986) and ana- 
lytical works (Bernardeau 1994). At sufficiently large radii 
where 5 is sufficiently small (S < 10) the fractional dispersion 
becomes small. At small radii, not only the dispersion be- 
comes larger, causing the mean value of S for a given Si to be 
somewhat different from S(Si) (expression (Q])), but also the 
effect of non-radial motions and velocity dispersions starts 
to dominate. So studying the dark matter profile where the 
uncertainty of the evolution is small, we may check that the 
initial conditions we use are correct to a high degree of accu- 
racy. With these conditions, as described by the most accurate 
probability distribution, P2(Si,q), we should be able to obtain 
very accurate predictions for all possible definitions (typical, 
mean...) of density and velocity profile. That is, there is no 
room for remodeling; all properties must be explained with 
one and the same P2(Si,q). Any residual discrepancy should 
be explained by triaxiality, non-radial motions and velocity 
dispersions, as we shall show in a future work, although, from 
the results shown in this work, the last two effects dominate at 
least up to roughly 5 virial radius. Going a step further, if we 
adequately take into account the neglected dynamical factors, 
the same initial conditions should explain the profile at any 
radii. In this way, we could be able to explain the dark matter 
profiles at least down to a virial radii understanding the role 
played by the initial conditions and by the different processes 
relevant to the evolution. 
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APPENDIX 
DERIVATION OF P 2 (S L , Q) 

The field of initial density fluctuations linearly extrapolated to the present Si{q) is a random Gaussian field. Any set of quantities 
obtained through a linear operation on a Gaussian field follows a Gaussian multivariate distribution. The quantities we are inter- 
ested in are the average of the linear field within three concentric spheres centered at a randomly chosen point with Lagrangian 
coordinate q. The Lagrangian radius of these spheres are Q, j (Q + q), q (note that q is not the norm of q) and we represent the 
average of the linear field within them by x\(q), x 2 (q), x^{q) respectively. These are three random Gaussian fields, but their values 
at a randomly chosen q we represent simply by x\, x 2 , X3. Their joint distribution is a Gaussian tri-variate: 

P(Xi,X 2 ,X 3 ) = - 

(2irp\detC\ 



jjXiXji (C)ij=<x iX j> 



' j j 



where C -1 is the inverse matrix of C, whose diagonal elements are the variances of the x's and the other elements are their 
correlations. 

Using the definitions of the probability distribution for x 2 conditioned to a given value of x\, P(x 2 \x\), and of the distribution 
for^3 conditioned to given values of x\, X2, P(xt,\x\ ,x 2 ), we may write: 

P(x u x 2 ,x 3 )=P(x l ) P(x 2 \xi) P(x 3 \x u x 2 ) (Al) 
Now, by grouping terms in \ in an appropriate manner we infer that: 



P(xi) = 



_ A i 
e ^ 

V27T <7\ 



e 2 « 

P(X 2 \ X1 ): 



2n gi 

: of (1-4); P = c n — 



P(x 3 \xi,x 2 )= 6 11' - (A2) 



2tt a' 

12— 2( 1 + 2 Cl2 C13 C 2 3-(c 2 l2 + C 2 n + cj 3 )) C l2 C 23 -C l3 cr 3 C12C13-C23 CT 3 

0=0% 3 ; A = -= ; a = -= 

(l-cf 2 ) l-cf 2 CT! l-cf 2 CT2 

<XiXj> 

c { j = ; ui = <j(qi) 

qi = Q; <72=^(2+<?); qi = q 

a(q) is the linear RMS density fluctuation as a function of Lagrangian radii (see expression dU). < x, Xj > is what we called 
Gij in expression ©, with qi, qj in the place of q, Q: 

1 1"°° 

<Jij = (Tij(qd = ^2 J \S k | 2 W T {qik) W T ( qj k) k 2 dk (A3) 
We may now obtain the distribution for X3 conditioned to X\ = 5 V! >, x 2 < 5 V i r : 

n , , , r . pfe 1 jcop(x 3 |x 1,^2)^2 

P(x 3 |jh = ^,„v,X2 < 6 vir ) = 00 (A4) 

P(x 2 \xi) dx 2 

Note that P(xi) cancels. 

Using ( IA2I ) and rearranging terms in the exponent we may write: 
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P(x 2 \x l )P(x 3 |xi,x 2 ) : 



1 



2ir g2 a' 



exp 



/ n \2 n2 1 (Px\.Bxi ABx\ \2 

(x^-AxxY P L x\ yTJ+TP? TP?> 



x exp 



2 a' 



f B 2 J_ 
\2a^ + Yg 



2g 



X2- 



2g + 2<t' 2 2cr' 2 ' 



2<r" 2g 



X2 appears only in the second factor and both integrals in (!A4t can directly be expressed in terms of the error complementary 
function: 

P(x-i\x x = S vir ,x 2 < 8 v ir) = — I : i — : exp 



\-\erfc[& 



2y/irg? L? a' 



(x 3 -AS rir ) 2 P 2 6 2 ir | U 2 
2g + L 



2 a 



12 



1 



1 erfc 

2 J 



S,ir(l-P) 



U ■ 



S vir P + Bx 3 ABS Vi 



2g 
L = 



2 a 1 
B 2 J_ 

2 a' 2 + 2g 



2 a 1 



Pi{5i,q), as defined in the text, is given by: 

Pi(Si,q) = P(x 3 = 6i\x\ = S vir ,x 2 < 5 vir ) 
So, P2(Si,q), which is simply P2{S/,q) with the restriction: 6/ < <5 VI > and the corresponding normalization, is given by: 



P(Xi=Si | xi=(5,, v ,j:2<(5 v ; r ) - f r r 
fE „, J - ] 7 TT~\ 11 °1 ^ °vir 

J '" P(-*3=°f X\=S vir ,X-,<S vlr ) 







if & >«y v 



(A5) 



Note that the right hand side in the case 8i > <5 n > depends on q through er,, c,y. The se quantities may be computed directly with 
arbitrary precision evaluating the integral entering their definition (expression (IA21 >). However, to be able to obtain efficiently 
P2(S, r), we need accurate fits to these quantities as a function of q. It is very difficult, however, to obtain consistent approximation 
for these quantities: very small errors in c (/ may lead to inconsistent result, for example, negative values for a' 2 . So, it is more 
expedient to fit directly the quantities {A,B,P,a' 2 ) where the er,-, cy enter. We find the following fits: 



A = 



-1.386 b I (p 2 -l : 



<73 
CTl ' 



B= 1 + 



(73 . 

C2 ' 



= -0.475 b(q?-l) 



P = e 



(A6) 



^-'2 _ 

a = a 3 



6.63 1(T 2 b 4 



((f) 2 " 1 ) 4 



where b is as defined in expression Q (b = 0.2544 for M = 3 x 10 12 /r'M ). This fit is valid in the range 3 x 10 11 to 
3 x 10 13 h~ l M Q . Ou tside this range it may be necessary using the full expressions (IA2t or a different fit to them. 
Using expression iA5i in expression d26l ) we find: 



P 2 (5,r) = -—G(S,q = Kl + 5)i) 



(A7) 



with: 



y DO 

G(6,q)= / P 2 (S l ,q = r(l 
JS,(S) 



+sy>)d8i 



where the dependence on q enters through the coefficients. 

Inserting expression ( IA7t into expressions (fl4] i and (l23l in the place of P(S,r), we may obtain the values of 5(r), V(r) (the 
average corresponding to S values between c* m ,„ and S max ) as given by our best approximation to the probability distribution for 5 
at a given r. After integrating by parts we find: 



S(r)-. 



G(<5, r) d5-5 m i n - 6 max G(S max , r) 



f{6) G(S,r) d5+f{5 min )-f(5 max ) G(8 max ,r) 



(A8) 



where f'(5) stands for the derivative of f(S) (defined in expression (fJTJ) with respect to 5. 



12 



BETANCORT-R1JO ET AL. 



APPENDIX B 

Let P(6i , q) be the probability distribution for the linear density fluctuation, Si, at a fixed Lagrangian distance, q, from the center 
of an object and let P(S, r) be the probability distribution for the actual density fluctuation, 8, at a fixed Eulerian distance, r, from 
the center of the same object. We shall show that the following relationship holds: 

[>00 [>O0 

/ P(6,r)dS= / P(S h q = r(l + 8o)hdSi (Bl) 
from which expressions (0 and d26| i immediately follow. 

If 5i takes the value 8i(8q) at Lagrangian radii r(l +80)3, it is clear that, by construction, S must take the value 60 at Eulerian 
radii, r. Now, if at Eulerian radii r, 5 takes a value , 5' , larger than Sq, its corresponding Lagrangian radii, q', must satisfy: 

q' = r(l + 5')i >r(l + <5 )^ =q 

and the value of Si at q', S' : = 5i{8'), must be larger than <5/(<5o)- But, if the linear profile is monotonically decreasing the value 
of the linear density fluctuation at q, 81, must be larger than the value at q', 8[, which is, ex hypothesis, larger than Si(Sq). So, 

whenever 8 is larger than 8q at r, 81 is larger than <5/(<5o) at q = r(l + Sq)t , hence, ( IBU follows. Note that, if P(8i,q) vanishes for 
81 > 8 V i r , the upper limit in the right hand side of (IBU may be set equal to <5,,,y. 
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